Van der Waals interaction in magnetic bilayer graphene nanoribbons 
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We study the interaction energy between two graphene nanoribbons by first principles calculations, 
including van der Waals interactions and spin polarization. For ultranarrow zigzag nanoribbons, the 
direct stacking is even more stable than the Bernal stacking, competing in energy for wider ribbons. 
This behavior is due to the magnetic interaction between edge states. We relate the reduction of 
the magnetization in zigzag nanoribbons with increasing ribbon width to the structural changes 
produced by the magnetic interaction, and we show that when deposited on a substrate, zigzag 
bilayer ribbons remain magnetic for larger widths. 

PACS numbers: 



I. INTRODUCTION 



Charge carriers in graphene follow a linear disper- 
sion relation close to the Fermi energy. For this reason, 
they are considered as massless fermions obeying Dirac's 
equation^ When several layers of graphene are piled up 
together, their electronic and transport properties can 
be dramatically modified, depending on the stacking ar- 
rangement and the number of layers^ There are sev- 
eral possible stacking arrangements in bilayer graphene, 
the most symmetric cases being direct (AA) and Bernal 
(AB) stackings. Most theoretical studies have focused 
on the AB stacking, because it is that of graphite, being 
the lowest energy configuration for the three-dimensional 
crystal^ However, the AA stacking has been observed in 
experiments on few-layer graphene, and it should also 
be considered in bilayer stackings4^— For example, the 
AA and AB stackings have been observed indistinctly at 
the graphene edges in samples grown on SiC4 In fact, 
bilayer graphene with AA stacking has also been synthe- 
sized, and observed by transmission electron microscopy 
(TEM) Due to the differences in the electronic proper- 
ties of bilayer AA and AB, the change between stackings 
by a relative displacement of the layers has even been 
proposed as the key mechanism for a switch device^r— 

Graphene nanoribbons are graphene strips of nanomet- 
ric width and arbitrary length, with electronic properties 
depending on their edges and widths! 10 ' 11 They are con- 
sidered as potential materials for future nanoelectron- 
ics because they can behave as metals or semiconduc- 
tors, making possible the design of electronic elements 
based solely on them. The simplest nanoribbon geome- 
tries are those with zigzag and armchair edges, which 
have been studied extensively,!! Other edges termina- 
tions are possible, but they can be mapped onto three 



basic types, the armchair being the only one without edge 
states . 12 ' 13 These localized states close to the Fermi en- 
ergy are responsible for the magnetic and transport prop- 
erties of zigzag graphene ribbons, and they are the origin 
of defect-related interface bands in graphene junctions ^ 
Within a simple tight-binding model, armchair graphene 
nanoribbons (AGNRs) can be either metallic or semi- 
conducting depending on their width^ whereas zigzag 
graphene nanoribbons (ZGNRs) are metallic with edge 
statesJ£ More realistic calculations yield all semicon- 
ducting armchair ribbons . 17 ' 18 With regard to ZGNRs, 
the inclusion of electronic interactions reveals a ferro- 
magnetic order of the magnetic moments at each edge, 
with an edge-edge antifcrromagnctic coupling that opens 
a small gap. In fact, this magnetic characteristic makes 
ZGNRs interesting for spintronic devicesi 19 ' 20 

In bilayer graphene nanoribbons (6-GNR) both, edges 
and stacking order, determine the electronic and mag- 
netic properties. Even though in few-layer samples there 
are multiple possibilities for the stacking arrangements, 
the majority of previous theoretical studies have focused 
on the AB stacking!^— The interaction between the 
edges of zigzag bilayer graphene ribbons determines the 
survival of magnetism. A combined first-principles and 
tight-binding approach was used to study the electronic 
properties in armchair and zigzag GNRi^ 4 - Because these 
authors find an important dependence on the functional 
employed, they fix their distance to graphite for zigzag 
GNR. Their magnetism is thus masked, as its survival 
depends on the layer-layer distance. An attempt to re- 
lax the edges was considered using a local spin density 
(LSD) approach within density functional theory^ They 
found that, for wide bilayer zigzag nanoribbons, the total 
magnetic moment is zero. 

Previous works do not consider van der Waals (vdW) 
forces. In order to relax bilayer ribbons, an explicit de- 
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scription of the vdW interaction must be included be- 
yond LSD. When these long-range interactions are in- 
cluded, the electronic densities between the layers are 
rearranged, and this yields variations on the interlayer 
distances. Such vdW interaction is included at a simple 
level in Ref. [25| and the edge magnetism disappears for 
small ribbon widths. However, we should note that vdW 
interaction is included in an effective way, modifying the 
atomic potentials. Other implementations using a fully 
nonlocal van der Waals density functional must thus be 
checked. 

In this work we study the properties of bilayer zigzag 
ribbons where all the edge carbon atoms are passivatcd 
by hydrogen. We include van der Waals dispersion forces 
with the fully non-local density functional recently pro- 
posed from first-principles within the family of func- 
tionals based on Ref. [28|, as recently factorized for 
efficiency^ In Section II we describe the computational 
details. We investigate the stability of Bernal and direct 
stackings in 6-GNR, focusing on the magnetic interaction 
between edges and on the interplay between magnetism 
and structural changes in narrow zigzag ribbons. Section 
III describes the systems studied, and shows our results, 
presenting the binding energies, and magnetic and struc- 
tural changes in zigzag bilayer nanoribbons. Our main 
conclusion is that direct stacking competes with Bernal 
stacking below a critical ribbon width, and we show that 
the magnetic coupling between edge states in the differ- 
ent ribbons plays a key role in such competition. Indeed, 
for ultranarrow ribbons, the direct stacking has the low- 
est total energy and largest binding energy. Furthermore, 
the structural distorsion at the edges due to this inter- 
action makes the magnetization negligible in bilayer rib- 
bons, causing metallization. However, when deposited 
on a substrate, the structural deformation is reduced, 
thus maintaining the edge magnetism for larger ribbon 
widths. We finish with a brief summary in Section IV. 



II. COMPUTATIONAL DETAILS 

First principles calculations are performed using the 
SIESTA code with spin polarization^ We use the van 
der Waals functional parametrized by Lee et al. (vdW- 
DF2) which is a second version of the original vdW-DF 
functional by Dion and coworkers^ The factorization 
proposed in Ref. represents a very substantial effi- 
ciency improvement in the evaluation of the exchange- 
correlation potential and energy, thus enabling first- 
principles van der Waals calculations for any system 
accessible to usual generalized gradient approximations 
GGAs. We check that the interlayer space in graphite is 
in agreement with previous calculations^ The results 
presented below have been performed using the func- 
tional vdW-DF2, but we have checked that other func- 
tionals implemented in the SIESTA code preserve the 
main features found employing vdW-DF2. Our choice of 
functional is motivated by the fact that vdW-DF2 gives 



more realistic binding energies for bilayer graphene when 
compared to experimental works. The electron- ion inter- 
actions use norm-conserving nonlocal Troullier-Martins 
pseudopotentials^I generated with the atomic configura- 
tion [He]2s 2 2p 2 taken as reference with a radius cutoff 
of 1.25 A for s, p, d and / orbitals. Spin polarized cal- 
culations normally require a fine sampling of the Bril- 
louin zone, which we performed with a Monkhorst-Pack 
scheme of 30 x 1 x I fc-points^ 2 - The real-space grid for 
matrix-element computations^ uses an energy cutoff of 
350 Ry. The structure was relaxed by conjugate gra- 
dients optimization until forces were smaller than 0.01 
eV/A. Periodic boundary conditions were applied, so we 
use large enough supercell parameters (15 A) in the di- 
rections perpendicular to the ribbon's long axis to avoid 
spurious interactions between adjacent ribbons. All the 
carbon atoms at the edges are passivated by hydrogen. 

III. RESULTS 

A. Ingredients: monolayer zigzag nanoribbons 

Before undertaking the calculation of bilayer nanorib- 
bons, we have first verified that our approach gives rea- 
sonable results for monolayer zigzag nanoribbons. For 
these, the key parameters to determine the electronic be- 
havior at the Fermi energy are both the edge shape and 
the ribbon width. We have used two initial magnetic 
configurations for the edges of ZGNR: either ferromag- 
netic {fm), with aligned spin polarizations, or antiferro- 
magnetic (afm), with antiparallel spin polarizations. No- 
tice that all atoms in the same edge are ferromagnetic 
coupled.— We have performed calculations of ZGNRs 
with several widths and both afm or fm ordcrings. We 
found that the afm order is always more stable than the 
fm, and that their energy difference decreases with in- 
creasing ribbon width, as in previous calculations.— 

B. Bilayer zigzag nanoribbons 

1. Binding energies and stable configurations 

As in infinite bilayer graphene, we have to look at dif- 
ferent stacking orders for bilayer graphene nanoribbons 
(&-ZGNRs). Fig. [1] (a) shows the three stackings con- 
sidered in this work. The top panel of the figure depicts 
an example of direct (AA) stacking. Two types of AB 
stacking have to be considered, according to the relative 
position of their edges. The medium and bottom panels 
of FigQ] (a) show the so-called AB^ and AB Q stackings 
for zigzag ribbons. We identify the ribbon width by N, 
being the number of zigzag chains from edge to edgeJ^ 

As the edges of 6-ZGNRs have magnetization, we have 
to study both the intralayer and interlayer (i.e., layer- 
to-laycr) magnetic couplings in these ribbons. We study 
four possible magnetic configurations for all the stackings 
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FIG. 1: (Color online) (a) Schematic packing of zigzag bi- 
layer graphene nanoribbons: A A (above), ARp (middle), and 
AB Q (below) . The dark gray layer corresponds to the bottom 
layer and the light gray (blue) layer correspond to the up- 
per layer. Hydrogen atoms are denoted by white balls. The 
ribbon width N is given by the number of dimers (zigzag 
chains) from edge to edge, (b) Magnetic configurations of 
the edge states in 6-ZGNR. Interlayer (intralayer) coupling 
can be either ferromagnetic [FM(/m)] or antiferromagnetic 
[AFM(flfm)]. Notice that one single edge is always ferromag- 
netic, i.e., all the spins along the same edge are parallel. 



considered, as depicted in Fig. [T](b). In order to distin- 
guish in our notation between intralayer and interlayer 
coupling, we use capital letters for the layer-to-layer cou- 
pling, and lower-case letters to label the intralayer cou- 
pling. The AFM- afm configuration (upper left diagram 
of Fig. [l](b)) thus has both antiferromagnetic intralayer 
and interlayer coupling. The bottom-left diagram, the 
FM- afm, shows two a/m-coupled GNR layers with FM 
interlayer coupling. The third and fourth configurations, 
AFM-/m and FM-fm, shown in the top and bottom right 
diagrams of Fig. [T] (b) , have both fm intralayer coupling 



with AFM or FM interlayer coupling respectively. 

We start from one of the four initial spin configura- 
tions for 6-ZGNRs described above in a range of dis- 
tances, and then we relax the minimum geometry so that 
the system evolves in principle towards similar magnetic 
configurations with lower energy. The converged solu- 
tions for different initial guesses are very close in energy. 
Each converged magnetic configuration can be viewed as 
a possible metastable solution.— It is likely that exter- 
nal conditions, such as magnetic and electric fields, can 
stabilize the system into in a configuration different from 
the energy minimum. 

However, for the AA stacking, the FM-(afm, fm) mag- 
netic initial guesses do not yield a stable solution: as 
the layers become closer, the electron density flips dur- 
ing self-consistency to the AFM ground state. The same 
happens for the AB stacking, where the FM- afm cases 
flip to AFM-afm solutions. Since for the AB stackings 
the atoms of an edge are not exactly on top of the atoms 
of the other, we obtained a large number of metastable 
magnetic alignments. 

From the total energies, we calculate the binding en- 
ergy (BE) as the difference between the energy of the 
coupled bilayer and the two isolated monolayers in the 
most stable configuration, i.e., the antiferromagnetic 
(afm) alignment^ The binding energy is related to the 
strength of the layer-layer interaction; it is given in meV 
per atom in a laycr j 35 ' 36 Figure [2] shows the binding en- 
ergy of &-ZGNRs with widths N ranging from 2 to 10, 
for several magnetic configurations between edges. The 
binding energies that we obtain for bilayer graphene are 
43.3 meV/atom and 50.6 meV/atom for AA and AB 
stacking, respectively. These values are in good agree- 
ment with the ones obtained from experimental works, 
namely, the exfoliation of graphite determines a binding 
energy between graphene layers of about 43 meV/atom^ 
and that estimated for the separation of polyaromatic hy- 
drocarbons is about 52 meV/atom.— As expected, we see 
in Fig. [2] that the binding energy increases in absolute 
value with the nanoribbon width. Interestingly, the in- 
crease is not monotonous and for certain widths it shows 
an enhanced stability, see for instance the N = 4 case. 

We find that the largest binding energy and most sta- 
ble configuration for > 4 is the AB Q in the AFM- 
afm configuration. Remarkably, for ultranarrow ribbons, 
N = 2,4, the AA stacking with AFM-afm coupling is 
more favorable, albeit with a very close value of the BE 
to that of the AB Q in the AFM- afm magnetic ordering. 

We now analyze in more detail the role of magnetic 
configurations on the stacking of ribbons. The cases AB a 
and AB^j which have fm intralayer coupling are in the 
same energetic range, with the binding energy of AB^ 
larger than AB a for all widths. With the exception of the 
ultranarrow widths, all the AB Q cases, as well as the AB^ 
with either FM or fm couplings, are rather close as to 
the binding energies. This shows a relationship between 
stackings and magnetic configurations of the edges. 

To elucidate the role of the magnetic interactions be- 
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FIG. 2: (Color online) Binding energy (BE) of b-ZGNR as 
a function of the ribbon width N after full relaxation. It is 
noteworthy that the AB a and AA stackings have the largest 
absolute values of the binding energies. 



tween edges, we have calculated the binding energies of 
ZGNRs on graphene. In such a case, because we are sup- 
pressing a ribbon, we are focusing on the edge-graphene 
interaction instead of the edge-edge interlayer interac- 
tion, and we only distinguish between AB and AA stack- 
ings. We find that the binding energy is lower in Bernal 
stacking for all the studied widths, as it has been previ- 
ously found in bulk graphite, where the Bernal stacking 
is more stable than the AA. This indicates that our find- 
ing on the greater stability and stronger binding energies 
for ultranarrow 4-ZGNR with AA stacking is related to 
the edge-edge interlayer coupling. 



2. Structural and magnetic changes: quenching of the edge 
magnetic moments 

The differences in binding energies between stack- 
ing orderings can be related to structural and magnetic 
changes in the bilayer ribbons. In our fully relaxed sim- 
ulations for 6-ZGNR, all the AB^ cases, and most of the 
AB a cases, the layers remain flat. Only the AA and 
AB a stackings with AFM-a/ra coupling change their ge- 
ometrical structure. Fig. [3] (a) shows two examples for 
N = 10. Notice that the edges are bent inwards; in the 
AA case, the layers bend symmetrically, becoming con- 
vex at their center, and in the AB Q case, for which the 
ribbons are laterally displaced, the edges approach main- 
taining a flat central region. The converged geometries 
of the ribbons are distorted, but still they have relevant 
symmetries, which are preserved within tolerance (~ 0.02 
A): the AA stacking with AFM-a/m configuration shows 
a mirror symmetry and Ci rotation with an axis parallel 
to the ribbons, while AB Q with AFM-a/m shows only C2 



symmetry. 
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FIG. 3: (Color online) (a) Structural distortions of the zigzag 
ribbon with N — 10 for the most stable stackings AA and 
AB a in the AFM-a/m configurations after full relaxation. The 
systems with higher binding energy (in absolute value) consist 
of non-planar, distorted, graphene ribbons, (b) Interlayer dis- 
tances at the center h c (empty diamonds) and edges h e (full 
diamonds) for 6-ZGNRs with AA (black) and AB a [blue, gray] 
stackings in the AFM-a/m configuration as a function of the 
ribbon width. 

To quantify these distortions, in Fig. GD(b) we plot the 
distances between the two ribbons at their center, h c , 
and at their edges, h e , as a function of the bilayer width. 
The distances at the central part h c remain constant for 
the AB a series, whereas they show larger changes for the 
AA stackings. For very small widths, up to N = 4, the 
central distances h c for AA stackings are lower than for 
the AB cases. This is not what happens in bulk bilayer 
graphene, where the layer-to-layer distance is smaller for 
Bernal stacking than for AA. This is an indication of the 
strong interaction in these ultranarrow ribbons with AA 
stacking. For bilayer ribbons with N > 6, the behavior 
is as expected, i.e., with central interlayer distances h c 
smaller in the AB^ ribbons than for the AA cases. On 
the contrary, the distance between edges, h e , is notably 
different from h c when N > 6. This deviation can be as 
large as 0.6 A, indicating a strong edge-edge interaction. 
Note that for the AB a cases, these structural distortions 
are accompanied by a lateral sliding, but the values we 
find, e,g., 0.1 A for N = 10, arc much smaller than those 
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FIG. 4: (Color online) Magnetization of the 6-ZGNR as a 
function of the ribbon width N after full relaxation. 

reported previously^ These interlayer distances h c and 
h e show that for ultranarrow ribbons, up to N = 6, the 
bilayer ribbons behave rigidly, becoming more flexible for 
larger widths. 

We now focus on the changes of magnetization at the 
edges M due to the edge-edge interaction, as also shown 
in Fig. fJJ It is defined as M — N up — Ndown, where 
N up (Ndown) is the number of electrons with spin up 
(down) per edge atom. The total magnetization by 95% 
corresponds to p z orbitals. We find that the magnetic 
moments are mainly located at the edges and decay ex- 
ponentially when moving into the central part of the rib- 
bon, in agreement with previous resultsi 10 ' 39 Notice that 
the interlayer interaction between edges suppresses the 
site magnetization: for the planar cases, the magnetiza- 
tion value is about 0.25 fib- This is the case for all the 
magnetic configurations with AB^ stacking, and the AB Q 
ones with FM intralayer coupling, see Fig. |U However, 
for the AA and AB a stackings with AFM-a/m couplings, 
when the interlayer edge distance h e is reduced, a strong 
interaction appears between the p z orbitals at opposite 
edges, and the spin cloud evolves towards nonmagnetic 
configurations. In fact, the spin polarization for A > 10 
is almost quenched. 



3. Implications of the magnetic quenching for calculations 
and experiments 

Narrowing of gaps at large widths. These magnetic and 
structural changes are associated with other variations of 
the electronic properties of ribbons. The band structures 
of the 6-ZGNRs of width N = 8 for AA and AB a stack- 
ings in the AFM-a/m configuration are shown in Fig. [5] 
(a). The gaps of the ribbons with A A stacking are smaller 
than those of the AB a ones. When increasing the ribbon 
width N the gap narrows as ~ 1/N^ Note the sharp 
drop in the energy gap after N = 6, related to the sudden 
decrease of the edge magnetization and the subsequent 
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FIG. 5: (Color online) (a) Band structure of a 6-ZGNR with 
N = 8 in the most stable stackings, AA and AB Q , and AFM- 
a/m magnetic configuration. The Fermi energy is set to zero. 
Note the narrowing of the gap. (b) Gaps versus ribbon widths 
for the two most stable stackings and magnetic configurations, 
namely, A A and AB Q all with AFM- afm. 

metallization of the bilayer nanoribbon, with E g < 0.05 
eV for N = 8 in the AA stacking, barely visible in Fig. 
[5] (a) i 10 i 39 It should be noted that the nanoribbon gaps 
may be underestimated in a GGA calculation^ The use 
of other functionals, such as the hybrid Heyd-Scuseria- 
Ernzerhof^i would correct this effect; in any case, the gap 
decrease with increasing width and the abrupt jumps as- 
sociated with the magnetic changes will certainly hold 
in calculations employing other functionals, but with the 
quenching of magnetic moments taking place at larger 
widths. 

Magnetoelastic switching. Our results show a strong re- 
lationship between structure deformation and magnetic 
configuration. However, a question that remains is on 
the reversibility of structural changes with respect to the 
magnetic configuration, which can be relevant for exper- 
iments. To address it, as well as to corroborate the inter- 
play between magnetism and structural changes, we have 
chosen the case of N = 10 in the ground-state, i.e., stack- 
ing AB Q and AFM-a/m configuration. As it is shown in 
the previous Section, the ribbons in this structure are 
strongly curved. The application of a magnetic field per- 
pendicular to the layer flips the magnetic moments at the 
edges from an AFM to an FM interlayer configuration; 
when we flip the magnetic moments from AFM to FM 
interlayer coupling in the curved structure and we relax 
it, it converges to a planar geometry in the FM configu- 
ration. We consider this finding to be an indication of a 
change from curved to planar geometry driven by mag- 
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nctic fields. In principle, in a magneto-mechanical device 
based on these ribbons, one could control the edge defor- 
mation with magnetic fields, which in turn can produce 
other electronic changes such as gap narrowing. 




x(A) 

FIG. 6: (Color online) Geometry of the 6-ZGNR with N = 10 
deposited on graphene. We have chosen a cut perpendicular 
to the graphene plane in order to highlight the edge defor- 
mations. The local magnetic moments on the edge atoms are 
also depicted. Note that the magnetic moments are almost as 
large as in a single ribbon. 

Effect of a substrate. We have studied the magnetic in- 
teraction of the edges in bilayer ribbons when deposited 
on a graphene substrate. Figure [6] shows the geometric 
and magnetic structure of the 6-ZGNR with N = 10 de- 
posited on graphene in the AFM-a/m magnetic configu- 
ration with AB a stacking. The interlayer distance in the 
center of the structure is overestimated, as it also occurs 
in bulk graphite^; this is a consequence of the vdW-DF2 
functional used. The top nanoribbon is deformed most, 
while the lower nanoribbon is nearly planar, due to the 
competing interaction between the graphene layer and 
the top ribbon. Due to this flat geometry, the associ- 
ated magnetization at the edges has higher values than 
the those obtained for bilayer nanoribbons, close to the 
ones of a single strip. As the structural deformation of 
the sandwiched layer is impeded by the substrate inter- 
action, the magnetic quenching is also precluded. There- 
fore, bilayer nanoribbons on substrates will remain mag- 
netic for larger values of N than when suspended. Bilayer 
nanoribbons with widths about 20 nm can therefore act 
as a spintronic device, maybe not in the free standing 
geometry, but certainly on substrates. 



IV. SUMMARY 

Bilayer zigzag graphene nanoribbons have been stud- 
ied by first principles DFT calculations including a vdW- 
DF2 van der Waals functional. Four possible magnetic 
configurations have been explored for the three more 
symmetric stackings, the AA and two Bernal (AB Q and 
AB P ). 

Our results show that the AA stacking is more favor- 
able for ultranarrow ribbons in the AFM-a/m configura- 
tion, competing in energy with the Bernal AB Q for larger 
6-ZGNRs. The edge interaction bends their structure in- 
wards for the AA and AB a stackings with an intralayer 
and interlayer antiferromagnetic configuration, but this 
bending is reduced for the smallest widths. With increas- 
ing ribbon width, the structural deformation at the edge 
is larger, leading to a reduction of the edge magnetic mo- 
ments and the metallization of the 6-ZGNRs. A magnetic 
external field can modify the structural changes, flatten- 
ing the ribbons. We have also studied the effect of a 
graphene substrate. In this case Bernal stacking is more 
favorable, and the bilayer ribbons maintain their magne- 
tization for larger widths. This is due to the reduction 
of the structural deformation because of the graphene 
substrate. 



V. ACKNOWLEDGMENTS 

This work has been partially supported by the Spanish 
DGES under Grants No. FIS2009-08744 and FIS2010- 
19609-C02-02, the Basque Departamento de Education 
and the UPV/EHU (Grant No. IT-366-07), and the 
Nanoiker project (Grant No. IE1 1-304) under the 
ETORTEK program funded by the Basque Research 
Dcpartamcnt of Industry. L. C. and H. S. grate- 
fully acknowledge helpful discussions with Juan Manuel 
Rodriguez Puerta and the hospitality of the Donostia 
International Physics Center. Calculations were partly 
done using the Camgrid high-throughput facility of the 
University of Cambridge. H. S. acknowledges the hospi- 
tality of the Department of Earth Sciences of the Univer- 
sity of Cambridge. 



* Electronic address: hernan.santos@icmm.csic.es; Present 
address: Instituto de Ciencia de Materiales de Madrid, 
CSIC. 

1 A. H. Castro Neto, F. Guinea, N. M. R. Peres, K. S. 
Novoselov, and A. K. Geim, Rev. Mod. Phys. 81, 109 
(2009). 

2 K. F. Mak, J. Shan, and Tony F. Heinz, Phys. Rev. Lett. 
104, 176404 (2010). 



3 J.-C. Charlier, X. Gonze, and J.-P. Michenaud, Europhys. 
Lett., 28, 403 (1994). 

4 W. Norimatsu and M. Kusunoki, Phys. Rev. B 81, 
161410(R) (2010). 

5 J.-K. Lee, S.-C Lee, J.-P. Ahn, S.-C. Kim, J. I. B. Wilson, 
and P. John, J. Chem. Phys. 129, 234709 (2008). 

6 Z. Liu, K. Suenaga, P. J. F. Harris, and S. Iijima, Phys. 
Rev. Lett. 102, 015501 (2009). 



7 



7 T. Ohta, A. Bostwick, T. Seyller, K. Horn, and R. Roten- 
berg, Science 313, 951 (2006). 

8 J. W. Gonzalez, H. Santos, M. Pacheco, L. Chico, and L. 
Brey, Phys. Rev. B 81, 195406 (2010). 

9 J. W. Gonzalez, H. Santos, E. Prada, L. Brey, and L. 
Chico, Phys. Rev. B 83, 205402 (2011). 

10 Y. W. Son, M. L. Cohen, and S. G. Louie, Nature (London) 
444, 347 (2006). 

11 L. Brey, and H. A. Fertig, Phys. Rev. B 73, 235411 (2006). 

12 A. R. Akhmerov and C. W. J. Beenakker, Phys. Rev. B 
77, 085423 (2008). 

13 W. Jaskolski, A. Ayuela, M. Pelc, H. Santos, and L. Chico, 
Phys. Rev. B 83, 235424 (2011). 

14 H. Santos, A. Ayuela, W. Jaskolski, M. Pelc, and L. Chico, 
Phys. Rev. B 80, 035436 (2009). 

15 K. Nakada, M. Fujita, G. Dresselhaus, and M.S. Dressel- 
haus, Phys. Rev. B 54, 17954 (1996). 

16 Y.-W. Son, M. L. Cohen, and S. G. Louie, Phys. Rev. Lett. 
97, 216803 (2006). 

17 M. Ezawa, Phys. Rev. B 73, 045432 (2006). 

18 V. Barone, O. Hod and G. E. Scuseria, Nano Lett. 6, 2748 
(2006). 

19 F. Munoz-Rojas, J. Fernandez-Rossier, and J. J. Palacios, 
Phys. Rev. Lett. 102, 136810 (2009). 

20 H. Santos, L. Chico, and L. Brey, Phys. Rev. Lett. 103, 
086801 (2009). 

21 E. McCann and V. I. Fal'ko, Phys Rev. Lett. 96, 086805 
(2006). 

22 E. McCann, Phys. Rev. B 74, 161403(R) (2006). 

23 J. Nilsson, A. H. Castro Neto, N. M. R. Peres, and F. 
Guinea, Phys. Rev. B 73, 214418 (2006). 

24 B. Sahu, H. Min, A. H. MacDonald, and S. K. Banerjee, 
Phys. Rev. B, 78, 045404 (2008). 

25 M. P. Lima, A. Fazzio, and A. J. R. da Silva, Phys. Rev. 
B, 79, 153401 (2009). 

26 G. Kim, and S.-H. Jhi, Appl. Phys. Lett. 97, 263114 



(2010). 

27 K. Lee, E. D. Murray, L. Kong, B. I. Lundqvist, and D. C. 
Langreth, Phys. Rev. B. 82, 081101(R) (2010). 

28 M. Dion, H. Rydberg, E. Schroder, D. C. Langreth, and 
B. I. Lundqvist, Phys. Rev. Lett. 92, 246401 (2004). 

29 G. Roman-Perez and J. M. Soler, Phys. Rev. Lett. 103, 
096102 (2009). 

30 J. M. Soler, E. Artacho, J. D. Gale, A. Garcia, J. Junquera, 
P. Ordejon, and D. Sanchez-Portal, J. Phys.: Cond. Matt. 
14, 2745 (2002). 

31 N. Troullier and J. L. Martins, Phys. Rev. B 43, 1993 
(1991). 

32 H. J. Monkhorst and J. D. Park, Phys. Rev. B 13, 5188 
(1976). 

33 N. Marom, J. Bernstein, J. Garel, A. Tkatchenko, E. Jose- 
levich, L. Kronik, and O. Hod, Phys. Rev. Lett. 105, 
046801 (2010). 

34 The binding energy is defined as -Ebe = Stotal — 2_Eml, 
where Stotai is the total energy of the 6-GNR and -Eml is 
the total energy of an individual layer in the most stable 
magnetic configuration, i. e., afm. 

35 S. Grimme, C. Miick-Lichtenfeld, and J. Antony, J. Phys. 
Chem. C 111, 11199 (2007). 

36 R. Podeszwa, J. Chem. Phys. 132, 044704 (2010). 

37 L. A. Girifalco and R. A. Lad, J. Chem. Phys. 25, 693 
(1956). 

38 R. Zacharia, H. Ulbricht, and T. Hertel, Phys. Rev. B 69, 
155406 (2004). 

39 A. Mananes, F. Duque, A. Ayuela, M. J. Lopez, and J. A. 
Alonso, Phys. Rev. B 78, 035432 (2008). 

40 J. Heyd and G. E. Scuseria, J. Chem. Phys. 121, 1187 
(2004). 

41 J. Heyd, G. E. Scuseria and M. Ernzerhof. J. Chem. Phys. 
118, 8207 (2003). 



